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Abstract. We investigate analytically the propagation of 
linear waves in a three-dimensional, nonmagnetic, isother- 
mal atmosphere stratified in plane-parallel layers. The mo- 
tivation is to study oscillations in the non-magnetic chro- 
mosphere and to assess the limitations of one-dimensional 
simulations of the K2 V bright point phenomenon. 

We consider an impulsively excited acoustic distur- 
bance, emanating from a point source, and propagating 
outward as a spherical acoustic wave accompanied by an 
internal gravity wave. The waves amplify exponentially 
in the upward direction. A significant wave amplitude is 
therefore found only in a relatively narrow cone about the 
vertical. The amplitude of the wave decreases with time. 
Because of the lateral spread, the wave amplitude decays 
faster in 2D and 3D simulations than in ID. The initial 
pulse, which travels at the sound speed, carries most of the 
energy injected into the medium. Subsequent wave crests 
leave the source region at ever-increasing phase speed, but 
slow to the sound speed as they approach the head of the 
wave. 

Important conclusions from the 3D solution that were 
not anticipated from the plane-wave solution are: 

1. The bulk of the energy is emitted in the upward (and 
downward) direction; much less goes into the horizontal 
direction. 

2. The wave profile narrows from the initial pulse through 
the amplitude maxima in the wake of the pulse. 

As a consequence of both points, the shock-heated regions 
in the wake of the initial pulse would weaken in strength 
and shrink in size. 

3. The height at which a given wave amplitude is reached 
spreads outward from the symmetry axis of the distur- 
bance as the wave propagates upward. Thus the diameter 
of the shock-heated region would increase as the acoustic 
wave travels upward in the atmosphere. 
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1. Introduction 

The solar chromosphere exhibits two signatures setting it 
apart from the underlying photosphere: an emission spec- 
trum that characterizes its thermodynamic state, and os- 
cillations that characterize its dynamical state. Both sig- 
natures are present in oscillations in H2 V and K2v bright 
points, which arise in the emission peaks on the blue side 
of the central absorption features in the resonance lines 
of Ca II, the H and K lines in Fraunhofer's nomenclature; 
these features originate at a height of approximately 1 Mm 
above r = 1 in the nonmagnetic chromosphere. 

The H and K lines are the strongest lines in the visible 
chromospheric spectrum and they show the three-minute 
oscillations. For these reasons they have often been the 
objects of observational studies (e.g., Liu 1974; Cram & 
Dame 1983; Lites, Rutten & Kalkofen 1993), and theoret- 
ical investigations (e.g., Fleck & Schmitz 1991; Kalkofen 
et. al 1994, hereafter KRBM; Sutmann & Umschneider 
1995; Sutmann, Musielak & Ulmschneider 1998) with the 
aim of elucidating the excitation mechanism of the waves 
as well as the properties of the chromosphere. 

An empirical simulation by Carlsson & Stein (1997) 
employing a sophisticated radiation-hydrodynamic treat- 
ment showed that the characteristic features of the spec- 
trum of the H line emerging from the chromosphere could 
be predicted from the observed velocity spectrum in the 
photosphere on the basis of propagating acoustic waves. 
Their model took the photospheric velocity field observed 
by Lites et al. (1993) and compared the computed time- 
dependent emergent H line profile with the observed pro- 
file. While the model gave an emergent H line intensity 
during the cool phase of the wave, as well as an overall 



2 



G. Bodo et al.: Acoustic Waves in a Stratified Atmosphere 



temperature structure, that are strongly contradicted by 
the observations (Kalkofen, Ulmschneider & Avrett 1999), 
the intricate velocity and intensity variations in the line 
core during the bright phase of the wave are reproduced to 
high fidelity; only in the timing of the computed H2 V in- 
tensity variation relative to the photospheric motion (and 
probably also in the absolute intensity) are there larger 
discrepancies from the observations. The otherwise close 
agreement between simulations and observations is sur- 
prising since the simulations, like all other previous the- 
oretical studies, assume plane acoustic waves. The reality 
of the wave propagation must be more complicated, how- 
ever, suggesting an idealization of cylindrical symmetry. 
It is therefore interesting to investigate the dependence of 
the results on the geometry of the waves and the limita- 
tions of the assumed plane symmetry. 

We model the excitation of linear hydrodynamic waves 
and their outward propagation in a three-dimensional, 
stratified, isothermal atmosphere, assuming that the in- 
teraction that gives rise to an outward traveling wave can 
be modeled as a pressure pulse at some reference level in 
the photosphere. 

Unlike the ID case, which allowed only acoustic waves 
to propagate, the 3D problem admits both acoustic and in- 
ternal gravity waves. Our treatment of the hydrodynamic 
equations allows the separation of the two modes. We will 
focus mainly on the acoustic mode. 

Questions to be examined include (1) the lateral spread- 
ing of the energy as the wave propagates upward, (2) the 
decay of the energy flux from its high value in the initial 
pulse to later oscillations in the wake of the pulse and (3) 
the increase of the phase velocity from the sound speed 
at the head of the wave to infinite phase speeds in the 
asymptotic limit. 

The paper is structured as follows: We present the ba- 
sic hydrodynamic equations in Section 2, describe the nu- 
merical results from the Fourier solution in Section 3, give 
the asymptotic analysis for late times in Section 4 and 
summarize the findings in Section 5. The appendices con- 
tain equations used in the Fourier solution of the problem. 



2. The Basic Equations 

The linearized hydrodynamic equations for adiabatic fluc- 
tuations, expressing the conservation of mass, momentum 
and internal energy for a gravitationally stratified, isother- 
mal atmosphere can be written in the form: 



Po- 



dpi 
dt 



-V ■ {p vi) 
-Vpi - pig 



(1) 



dpi „ of dpi _ 

-j^ + vi ■ Vp = a 2 I + vx ■ Vp 



where pi, pi and Vi are, respectively, the perturbed gas 
pressure and density and the velocity; the subscript de- 
notes the unperturbed variables; g = {0, 0, g} is the grav- 
itational acceleration, a = VTPoTpo is the sound speed 
and 7 is the ratio of specific heats. The equilibrium values 
of unperturbed variables have the form: 

Po = Poo exp (-z/H) , po = poo exp (-z/H) ,v = 0, and 



const . 



where poo and poo are, respectively, the values of pressure 
and density at the reference height z = 0; and H = a 2 j^g 
is the pressure scale height. 

In order to simplify the calculations it is convenient to 
assume the height dependence (see Lamb 1932) 

^loce^, (2) 

and, as a consequence, 

pi , pi cx e- z ' 2H . (3) 

With these assumptions and measuring space in units of 
twice the pressure scale height, 2H, and time in units of 
uj^ 1 — P ac /27r, where P ac = 4nH/a is the cutoff period, 



equations <W become 



<9V 



— - - (A - 1)— - - 4- — 

dt* 1 ' dt 2 7 2 



A_l<^ = 



(4) 



where tp is any of the variables (||) or (||) in non-dimensional 
form, i.e., v — Vi/a, p = pi/po and p = pi/po, i.e. the 
'reduced' perturbations, and A and Aj_ are Laplacian op- 
erators. Note that the reduced quantities p and p bear a 
vertical amplification as exp (z/2H) even though the ac- 
tual perturbations pi and pi decay with z. Note also that 
Eq. (|J) reduces to the Klein-Gordon equation in the one- 
dimensional limit (see, e.g., KRBM). 

We perform a plane-wave analysis assuming that all 
variables are proportional to exp[i(fc-r — tot)] . Equation (Q) 
then yields the dispersion relation 

,7-1 



^-(l + k 2 ) u/+4- 
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, ki = 0, 



(5) 



(Bray & Loughhead 1974), where k 2 = k 2 + k\ is in units 
of 1/2H and u> in units of w ac . Eq. (||) can be solved to 
obtain 
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where the plus sign is for the two roots for acoustic waves 
and the minus sign for the two roots for internal gravity 
waves. 

The solution of system ([!]), subject to given initial con- 
ditions at t = 0, can be written in the form of Fourier 
Integrals: 



p= / Aj(k) exp[i(k ■ r — ujjt)]dk (7) 

7=1,4"' 

/ A J -(fc)P(fe,w j )exp[i(fe-r-Wji)]dfc (8) 

1 A J 
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V I = 



^2 I Aj(k)G(k,u]j)exp[i(k-r -LJjtydk (9) 

=1,4'' 

/ AjityKfcu^expliik-r-WjtJldk (10) 

1 A J 



j=l,4 



The summation goes over four modes corresponding to the 
four roots of the dispersion relation (||), i.e., both acoustic 
and internal gravity modes travelling in two directions. 

In order to obtain the values of F, G, K of Eqs. [8^101, 
i.e., the relation between density, velocities and pressure in 
the plane- wave solutions, we have Fourier analyzed system 
(|l|), derived a set of algebraic equations (Appendix A), 
and determined the amplitudes Aj(k) by matching the 
given initial conditions at t = (Appendix B). As initial 
conditions we have considered a pressure pulse; the set of 
conditions becomes: 

p(r, t = 0) = Sp exp (— r 2 /r$) 

p(r,t = Q) = 

v z (r,t= 0) = 

k± -v± (r,t = 0) = 

where dp is the amplitude of the initial perturbation and 
ro is the spatial width of the Gaussian. Note that this 
four conditions completely determine the solution since 
we have fourth order equation. The analytical form for 
the amplitudes Aj(k) is given in Appendix B. 



3. Analysis of the Results and Asymptotic 
Solution 

Before we solve these equations it is instructive to con- 
sider briefly the known cases of waves in one-dimensional 
media. For a plane wave in a stratified atmosphere the 
hydrodynamic equations admit only acoustic waves, but 
no internal gravity waves. The solution for impulsive ex- 
citation of a disturbance (see KRBM) gives an upward 
traveling pulse that amplifies exponentially with height; 
it is followed by a wake that oscillates at the acoustic cut- 
off period. By contrast, a plane acoustic wave due to an 
impulse in a homogeneous medium, for which the atmo- 
sphere of our every-day experience is a natural approxi- 
mation, shows only the signal, but no wake. It is interest- 
ing to investigate the propagation of a disturbance in a 
3D medium and to highlight the differences with the ID 
cases. In order to gain some preliminary insight we study 
the asymptotic behavior of the solution in 3D. 



3.1. Asymptotic analysis 

Following Whitham (1974) we describe the solution in 
terms of a slowly varying wavetrain (i.e., with little vari- 
ation in a typical wavelength and period), writing, for ex- 
ample, the pressure perturbation as 
P (r,t)~A(r,t)exp{i6(r,t)), (11) 




2 
k 



, k 



Fig. 1. Group velocity vs vertical (solid line) and horizon- 
tal (dashed line) wavenumbers. 



where the amplitude A and the phase 9 are slowly varying 
functions of position and time. We can then define a local 
wavenumber, k = V0, and a local frequency, to = —d9/dt, 
to obtain for the wavenumber k the equation of motion 
dk 

— + (v g -V)k = , (12) 

(see Eq. 11.44 in Whitham, 1974) where v g is the corre- 
sponding group velocity vector. The equation is valid in 
the limit t — > oo. 

According to Eq. S3), the wavenumber k is constant 
along group lines, which are defined by the equation, and 
each value of k propagates with the corresponding con- 
stant group velocity v g (k). From the initial pulse, given 
by the superposition of modes encompassing the whole 
spectrum of wavenumbers, each wavenumber then propa- 
gates with its own, corresponding group velocity. A par- 
ticular wavenumber found at position (x, z) at time t can 
be obtained by solving the pair of equations 

x = [v s (k z , k±)]±t (13) 
z=[v g (k z ,k ± )} z t. (14) 

for the wavenumber components. The group velocity com- 
ponents in the vertical and horizontal directions are given 

by 



(Vg)z = 



(V 



duj 


-- h. 


dk z 




dui 


k± 


dk ± 


2oj 



1 



V(fc 2 + l) 2 -16^( 7 -l)/ 7 2 
fc 2 + l-8( 7 -l)/ 7 2 



(15) 



^(k 2 + l) 2 -16kl( 1 -l)/ 1 2 

where u is a root of Eq. (^|) with the plus sign. 

It is evident, from the dispersion relation (^) or from 
the system (^) - ([!(]), that the complete solution is given 
by the superposition of acoustic and gravity waves. Here 
we confine our analysis to the acoustic waves. 

The group velocity as a function of the moduli of the 
wavenumbers for the vertical and horizontal directions is 



(16) 
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Fig. 2. Vertical (upper panel) and horizontal (lower panel) 
distributions of wavenumbers. 

plotted in Fig. [ij which shows that the behavior of the 
two curves is different at small and intermediate wave- 
numbers: in the horizontal direction the group velocity 
remains small up to k ~ 0.7 and then increases rapidly to 
the asymptotic value of unity, and in the vertical direction 
the growth begins steeply at k = and then continues 
more gradually to the asymptotic value. 

The evolution of the distribution of wavenumbers in 
space is described by Eqs. (|13f ) and fll4| ), with the group 
velocity v s given by Eqs. (|l5|) and (|16|). These equations 
can be solved for the components k z and k± in the verti- 
cal and horizontal directions as functions of z/t and x/t; 
the solutions are shown in Fig. ||. The two panels reflect 
the different behavior, described above, of the group ve- 
locity as a function of wavenumber; in both cases, small 
wavenumbers, which have low group velocity, are found 
close to the origin, while large wavenumbers, which have 
group velocity close to the speed of sound, are found near 
the pulse. At intermediate distances, the distribution of 
k values is much broader in the vertical than in the hori- 
zontal direction. This has immediate consequences for the 
nature of the solution: since the value of the wavenum- 
ber found at a particular position implies also the wavele- 



length of the oscillation, one may expect (see Fig. g below) 
that in the horizontal direction the oscillations behind the 
front have almost constant wavelength, while in the verti- 
cal direction the wavelength gradually increases from the 
head of the wave towards the origin. 

Combining the phase velocity from the dispersion rela- 
tion with the distribution of wavenumbers obtained from 
the evolution equations (|l3| ) and (|lj) one can follow the 
propagation of points of constant phase, for example of a 
maximum. When a particular maximum is located close 
to the origin (x, z = 0) it is characterized by low val- 
ues of the wavenumber, and thus high phase velocity. As 
it moves away from the origin, its local wavenumber as 
well as its frequency increase (in the vertical direction, 
u} 2 — k 2 + 1), while its propagation velocity decreases 
(Fig. |^). Thus, as a particular maximum travels away from 
the origin and towards the head of the wave it is increas- 
ingly characterized by high-frequency components (which 
become important for shock formation in the nonlinear 
regime). Fig. ^ presents this behavior in more detail; it 
shows the phase velocity of three consecutive maxima as 
functions of distance from the source for the vertical (up- 
per panel) and horizontal (lower panel) directions. In the 
horizontal direction the phase velocity is practically equal 
to the sound speed nearly throughout the whole region, 
and for all maxima. In the vertical direction the behav- 
ior is more complicated: Except for the head of the wave 
(the first maximum in Fig. |[ upper panel), which prop- 
agates at the sound speed, the later maxima have high 
phase velocity over increasingly extended height ranges. 
They therefore travel at increased phase velocity (Fig. |[ 
upper panel; see also KRBM Fig. 2), especially near the 
origin (where i> p h — ► oo as z — > 0), until they approach the 
head of the wave. As a consequence of this difference in 
behavior between the two directions, the surfaces of con- 
stant phase behind the wave front travel a longer distance 
in the vertical than in the horizontal direction and thus 
acquire oval shape (see Fig. || below). 

3.2. Analysis of the results 

We calculate the Fourier integral (7) numerically, taking 
the symmetry with respect to the z-axis into account, and 
interpret the results with the aid of the asymptotic anal- 
ysis. Fig. ^ shows a representative solution for a pressure 
pulse at x=z=t=0, displaying the reduced pressure dis- 
tribution, i.e., without the vertical amplification factor 
exp (z/2H), in a vertical plane containing the origin, at a 
time t corresponding to about three times the acoustic cut- 
off period. As discussed above, the solution is given by the 
superposition of acoustic and internal gravity modes. The 
pressure variation due to internal gravity waves, which 
appears as radial stripes, remains confined to the vicinity 
of the origin since their group velocity is lower than that 
of acoustic waves. It is also highly anisotropic: internal 
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Fig. 3. Phase velocity for three consecutive maxima (solid: 
first maximum, dashed and dash-dotted: second and third 
maxima) in the vertical (upper panel) and horizontal 
(lower panel) directions. 




Fig. 4. Acoustic and internal gravity waves at t = 3 times 
the acoustic cutoff period, in a vertical plane through the 
origin. Wave crests and troughs are displayed by bright 
and dark shading, respectively. 



gravity waves are excluded from the purely vertical direc- 
tion, and in the horizontal direction their group velocity 
reaches only about 60% of the sound speed. Because of 
their higher group velocity acoustic waves become more 
prominent at greater distance from the origin. They are 
seen alone in Fig. g, which shows that the head of the 
wave forms a spherical surface, with radius equal to a x t. 

Recall that in the one-dimensional case the solution 
consists of a wave front followed by a wake that oscillates 
at a period approximately equal to the acoustic period 
(u> = 1). In the three-dimensional case this property is 
found again in the vertical direction (compare Fig. ^, up- 
per panel, and KRBM, Fig. 2), but now also in the hori- 
zontal direction, albeit with much reduced amplitude. 

The oval shape of the pressure extrema is related to 
the distribution of fc— values in space (Fig. |2[). As seen in 
Fig. ^ (lower panel), the pressure amplitude in the hori- 
zontal direction has wave crests that are nearly equidis- 
tant, consistent with the narrow distribution of A;— values 
and, consequently, wavelengths, through most of the range 
of x— values. By contrast, the distance between pressure 
extrema in the vertical direction increases with distance 
from the head of the wave due to the corresponding in- 
crease in the phase velocity, from the sound speed at the 
head of the wave to larger values near the origin (Fig. 2), 




Fig. 5. Acoustic waves at t — 3 times the acoustic cutoff 
period in the vertical plane. 
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Fig. 6. Vertical (upper panel) and horizontal (lower panel) 
cuts of the pressure perturbation amplitude (in arbitrary 
units) at t = 3 times the acoustic cutoff period. 



(and to infinitely large values in the asymptotic limit of 
t — > oo). 

Since the disturbance is due to a point source, the wave 
starts out as a spherical wave. The pulse retains that shape 
but as can be seen in Fig. [j| the amplitude of the pulse 
increases from the vertical direction (both positive and 
negative) towards the horizontal direction. 

The three-dimensional solution allows us to investigate 
the variation of the perturbation amplitude with direction 
relative to the vertical. We can look to this variation from 
two different perspectives: one possibility is that of con- 
sidering the variation of the amplitude of a maximum at 
a fixed time and the other is that of considering it at a 
fixed height (see Fig. ^|). In an astrophysical context this 
second case would correspond to looking at a fixed optical 
depth and therefore to a particular spectroscopic signature 
of the wave. In the first case (fixed time), the variation 
of the amplitude has two parts, the intrinsic variation of 
the reduced function, and the variation of the exponential 
factor, exp(z/2iJ), that compensates for the exponential 
dependence of the background mass density on height and 
converts the reduced variables into the physical variables. 



N 




10 15 20 25 
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Fig. 7. Upper panel represents the positions of a wave 
crest at different times and their intersections with height 
of 15H. Lower panel shows a wave crest whose apex has 
reached the height of 15H at time to. 



But at a fixed height the distinction between reduced and 
physical variables is irrelevant. 

We consider first the case of fixed height: Fig. || shows 
the variation of the amplitude of the pulse and of the first 
three maxima of the wave train as they cross the height 
of 15£f , as functions of the angle measured relative to the 
vertical direction. Note that, in the actual solar chromo- 
sphere, nonlinear effects would already be present at this 
height for typical initial perturbations; however, our mo- 
tivation is to study the basic linear behavior of the wave 
propagation. Since the various points of the pulse, for ex- 
ample, reach the target height of 15 H with position x and 
hence angle increasing with time (cf.. Fig. |^), the variation 
of the pressure amplitude reflects partly the spherical de- 
cay of the pulse with distance r from the origin (p ~ r -1 ), 
and partly the increase of the amplitude with angular dis- 
tance from the vertical, as seen in Fig. ^. The net effect is 
a relatively slow decrease of the pressure amplitude, which 
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Fig. 8. Angle dependence of the (reduced) pressure am- 
plitude relative to the vertical direction (0°). The values 
are for the fixed height of 15H and refer to the pulse ( — ), 

and the 1st ( ), 2nd (— • — ) and 3rd (—••■—) maxima 

of the wake. 



reaches a factor of 2 at 60°. The decay is much faster for 
the wave crests in the wake, reaching a factor of 10 at 60° 
for the first and second maxima, and at 40° for the third 
maximum. The wave profile evidently narrows consider- 
ably with time. 

A picture complementary to Fig. || is shown in Fig. [|, 
with the angle variation of the pressure perturbation at 
the times when the apices of the pulse and, subsequently, 
of the maxima in the wake reach the height of 15H . The 
upper panel gives the angle dependence of the reduced 
pressure (at the time when the apex reaches z = 15H), 
and the lower panel that of the physical pressure. The up- 
per panel confirms the impression from Fig. |^ that the (re- 
duced) pressure amplitude of the pulse grows with zenith 
angle. The crests of the wake weaken with separation from 
the vertical; they also have significantly lower amplitude. 
The influence of the exponential factor on the physical 
pressure is clearly seen in the lower panel, and it results 
in a much faster drop of the amplitudes with zenith an- 
gle. This is especially true for the maxima in the wake for 
which the oval shape causes a much faster decay than for 
the pulse and therefore a narrowing of the wave profile 
with the order of the maximum. For the pressure pulse 
(Fig. lower panel), the decay in a half-angle of 45° is 
by nearly two orders of magnitude, and for the maxima in 
the wake, by nearly three orders of magnitude. Note that 
this decay reflects to a large extent the shape of the wave 
crests and is much less dramatic when the points along the 
crests have reached the target height, as described above. 
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Fig. 9. Variation with zenith angle of the pulse 

the first three maxima of the wake ( ,— • — , — ■ • ■ — , 

resp.) at the instants when the apex for each feature has 
reached the target height of 15H; the upper panel shows 
the reduced pressure, the lower panel includes the expo- 
nential height factor. 

4. The energy propagation 

For the system of the linearized equations (Q) we can write 
the conservation equation for total energy in the form 



dE 



o. 



(17) 



where the total energy density E (in units of poa 2 ) is de 
fined as 

v ■ v H — - + | p 



E = 



1 



(18) 



7 2 2 7 - 1 
and the energy flux F as 

F=pv. (19) 

In expression (|l8|) for E the terms represent, respectively, 
kinetic energy, acoustic potential energy and gravitational 
potential energy. For the slowly varying wavetrain (11) 
we may consider energy and energy flux averaged over a 
period, denoted by £ and where the averaged energy 
flux may be written as 

T = v s £ (20) 

(Whitham 1974, Eq. 11.69). This (asymptotically correct) 
form of the energy flux shows the conservation of total 
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energy in any volume in ordinary space whose boundaries 
move with the group velocity according to Eqs. (|l3r|lT 
Eq. (|l^) can therefore be written as 
d£ 



Acoustic at x = 0, z = 2tt 



(21) 

One can also read Eq. noting that a positive diver- 
gence of the group lines yields a decay of the energy den- 
sity. Now, since £ oc \A\ 2 one can derive the temporal 
evolution equation for the amplitude \A\ 2 , 

(see Whitham 1974, Section 11.6), where n is the number 
of dimensions. Eq. (E^) yields 

\A\ oc t- n/2 . (23) 
Fig. [l(] shows the temporal decay of the reduced pres- 
sure amplitude of the wake after the pulse has reached the 
positions of z = 2tt on the vertical axis (top panel) and 
y = 2n on the horizontal axis (lower panel) . The top panel 
compares the decay law of i~ 3 / 2 with the actual behav- 
ior of the wake. The agreement is excellent even though 
the analytic result ( p3| ) is valid only for a slowly varying 
wavetrain. Note the (real) time delay of ~ P ac /2 in the 
horizontal direction with respect to the vertical, a conse- 
quence of the phase velocity of the wake that is smaller in 
the horizontal direction (see Fig. ||). Note also that the am- 
plitude of the maxima in the horizontal direction is lower 
by a factor of 2 to 3. This follows from the properties of 
the energy propagation, as will be discussed below. 

At t = the energy is spatially concentrated in the ini- 
tial pulse, and as time elapses, the energy is dispersed into 
a wavetrain and the wavenumber distribution is spread 
out in ordinary space; both occur with the group veloc- 
ity, as described by Eqs. ( |l2[) and (|2l|). Thus the energy 
in any volume in wavenumber space remains fixed as the 
wavenumbers spread in physical space. Consequently 
<S(k(x, y, z, t)) dk x dk y dk z = E(x,y, z,t) dxdydz , 
where E(x,y, z,t) is the energy per unit volume in ordi- 
nary space, estimated at the point (x,y,z) and time t, 
and (£(fc(x, y, z, t)) is the energy per unit volume in the 
corresponding wavenumber space, where k(x, y, z) is the 
wavenumber at position (x, y, z) and time t. Now, volumes 
transform as 



dk x dkydk z — 
where 
J 



d{x,y,z) 



d(k X7 ky, k z 



dxdydz , 



d{x 7 y,z) 



d(k x , ky, k z 

is the Jacobian of the tranformation between wavenum- 
bers and ordinary space; J depends on the initial dis- 
tribution and the expansion or contraction of volumes in 
the mapping from wavenumber space to ordinary space, 
defined by the system ([lj|-[l4|). The energy density can 
therefore be written as 

<£(k(x,y,z,t)) 



3 F 



1; 



a; 




■1 F 
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0) 
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(24) 



Fig. 10. Pressure amplitude of the wake (in arbitrary 
units) vs time at x = and z = 2ir (top panel) and 
at x = 2tt and z = 0. The pulse has been removed; the 
amplitude decay of £~ 3 / 2 is given as a comparison in the 
top panel. 



From the above discussion it is clear why the ampli- 
tude of the maxima following the pulse is much smaller in 
the horizontal than in the vertical direction. The behavior 
of the group velocity and the distribution of k— values in 
space (see Figs, [l] and ||) show that in the horizontal di- 
rection a narrow range of k values near k = 1 is dispersed 
over a large portion of physical space, while the dispersion 
is much lower in the vertical direction; energy is therefore 
much more diluted in the horizontal than in the vertical 
direction and, as a consequence, the maxima have much 
lower amplitude. 

This behavior becomes even more evident from the 
fractional energy contained in a cone of infinitesimal open- 
ing angle when its symmetry axis points in different direc- 
tions. Fig. ^ shows this fractional energy at the time when 
the front has reached the height of 15H , normalized to the 
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distance from the front 

Fig. 11. The fractional energy contained inside a cone vs 
distance from the front, for a cone whose axis is inclined 

at 0° (— ), 45° (• • •), and 90° ( ) with respect to the 

vertical. 

total energy contained in a vertical cone, against distance 
measured from the front of a wave. Note that in the hori- 
zontal direction the energy rises from zero at the front to 
the asymptotic value in ~ 0.5ff. This implies that most 
of the wave energy is contained in a thin layer, of ~ 0.5ff 
width, following the pulse. But in the vertical direction 
about 80% of the energy is contained within a layer with 
a thickness of ~ 5ff extending in height from lOff to 15H. 
The wake therefore contains a significant fraction of the 
energy. Approximately 30% of the energy is in the pulse 
itself, i.e., within a thin layer of ~ 0.5-ff at the head of the 
wave. Note also the sharp drop of the energy with polar 
angle in Fig. |ll|. Evidently the bulk of the energy goes into 
the vertical direction. 

5. Conclusions 

We have carried out a linear analysis of the three-dimen- 
sional wave propagation in a stratified, isothermal atmo- 
sphere. The motivation was to study the limitations of the 
generally assumed plane waves. Although the analysis was 
based on linearized equations and idealized medium our 
results may suggest the sense in which 3D and ID waves 
differ from one another in the nonlinear regime. 

A pulse generated in the 'photosphere' by a 'point 
source' propagates upward with exponential amplification. 
The source region has a Gaussian pressure perturbation 
with a full half-width of 0.4 x ff , where ff is the pressure 
scale height, for a diameter of 40 km. An initial perturba- 
tion of strength Sp/p ks 0.2 becomes nonlinear at a height 
of 10 x ff . So this strength of the source is near the min- 
imum needed for K2 V bright points. 

At 1 Mm above the source, which corresponds to the 
height of formation of the H and K lines of Ca II in the so- 
lar chromosphere, the pressure and velocity perturbations 



reach large values over a region with a horizontal extent of 
1 Mm. This is comparable to the size of Ca bright points. 

The initial pulse is followed by a wake in the upward di- 
rection at the acoustic cutoff period, approximately three 
minutes in the solar atmosphere. There is a wake also in 
the horizontal direction, with the same period but much 
lower amplitude. 

The energy of the wave is concentrated in the vertical 
direction: One quarter of the upward-propagating energy 
is contained within a cone with a half-angle of 30° about 
the vertical axis; that cone constitutes only about 13% of 
the volume of the hemisphere. The energy is concentrated 
also in a narrow layer behind the initial pulse: When the 
wave has reached a height of lOff above the source, 60% 
of the energy on the vertical axis is contained within the 
first two scale heights behind the apex. 

The height where a given magnitude of the pertur- 
bation is reached increases with the size of the region. 
Assuming that nonlinear conditions are reached on the 
vertical axis at a height of lOff, a 'bright point' with 
a diameter of 1 Mm requires that the wave at the edge 
travel another 1.5ff, which requires an additional 21 s (at 
a sound speed of 7 km/s), and for a diameter of 2 Mm, 
the additinal travel distance and time are 3ff and 45 s. 
Thus the 'bright point' grows from the center outward and 
upward. For the larger diameter and height, correspond- 
ingly higher layers in the chromosphere would form the 
spectrum. 

The amplitude of the oscillations behind the pulse weak- 
ens in strength and shrinks in size. As a consequence, the 
maximal intensities associated with bright points in the 
wake are weaker and smaller than those in the initial pulse. 

The maxima in the wake decay with time as v(t) oc 
d e p enc hng on the geometry. For a plane wave in a 
one-dimensional medium, n = 1, the velocity decays as 
i^ 1 / 2 ; it matches the plane wave solution of KRBM; for a 
line source in a two-dimensional medium, the decay is as 
i -1 ; and for a point source in 3D, it is as v{t) oc t~ 3 / 2 . The 
energy flux decays as v 2 {t). Thus the energy contained 
within the spherical volume of the wave decays in step 
with the increase in volume and hence the energy within 
is conserved. 

For comparison with observations of chromosphcric os- 
cillations it needs to be borne in mind that the initial 
atmosphere in our analytic solution is at rest. The numer- 
ical simulations by Carlsson & Stein (1997) compare well 
with observations only when the waves are launched into 
a disturbed atmosphere. 
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Appendix A: 



System (|l|) can be Fourier-analyzed in terms of Vi, pi, pi oc 
exp(±z/2//) exp[i(fe ■ r — wt)], with the ± sign according to 

(I. 

(Al) 
(A2) 



), to obtain a set of algebraic equations: 
— ik ■ vi = 

Pi 



iuipovi — ikpi + ( 



ibjpi — iua pi 



(g 



2H 
2H 



gpi)z = 



)poviz = , 



(A3) 



The linearized variables can be put in a non dimen sion al form, 
as described in Section 2, and Eqs. (Al), (A2) and (A3) can be 
solved to obtain p = Fp , v z = Gp , k±-v± — Kp, with F, G and 
K specifying relations between density, velocities and pressure 
in the plane-wave solutions as: 

j-i/ i \ k -\- 2ihz 1 / \ t \ 

F ( UJ > k ) = — —i 77TT ( A4 ) 



G(u,k) = — 



70; 



7a; 2 - 2(1 - ik z ) 
ijuj 2 k z - 2[kl + k 2 ) 



2(1 -ikz) 
K(u,k) 



fei 
701 



(A5) 
(A6) 
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Sutmann, G., Musielak, Z. E. & Ulmschneider, P. 1998, AA 
340, 556 

Whitham, G. B., 1974, Linear and nonlinear waves, John Wiley 
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Appendix B: 

The amplitudes Aj(k) for an initial pressure pulse are: 
k\ F{lo 2 ) 



Ai = Az = Sp exp (— 
A2 = A4 — —Sp exp ( 



kl'F(uo 2 ) -F(ui) 

k 2 F(lui) 

k 2 ' F(u 2 ) - F(loi) ' 



where fco = 2/ro, 0J1 is the solution of dispersion relation 
for acoustic modes and ui 2 is the solution for internal grav- 
ity modes. Ai and A3 are the amplitudes of acoustic modes 
while A 2 and Ai are the amplitude of internal gravity modes. 
For an initial velocity pulse we obtain: 

a a r / k 2 . K(ui 2 ) 
Ai = -A 3 = -Sp exp (--j) i— i 



k 2 ' K(uji)G(uj 2 ) - K(u2)G(ui) 



a a r 1 k 2 , K(uji) 

A 2 = -Ai = -Sp exp(- T j) 



k 2 ' K{uji)G{uj 2 ) - K{uj 2 )G(uji) 



The solution has then been obtained by performing numerically 
the Fourier integrals (p[|R|). 
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